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We confront, quantitatively, the theoretical description of the reaction-diffusion of a second order 
reaction to experiment. The reaction at work is Ca 2+ /CaGreen, and the reactor is a T-shaped 
microchannel, 10 fim deep, 200 /im wide, and 2 cm long. The experimental measurements are 
compared with the two-dimensional numerical simulation of the reaction-diffusion equations. We 
find good agreement between theory and experiment. From this study, one may propose a method 
of measurement of various quantities, such as the kinetic rate of the reaction, in conditions yet 
inaccessible to conventional methods. 



Chemical reaction processes taking place in spatially 
extended reactors with no active mixing are described 
by reaction-diffusion equations. In general, these equa- 
tions are non linear and they possess a rich variety of 
solutions, giving rise to many different spatio-temporal 
patterns. Reaction-diffusion equations have been known 
for more than a century but the first elementary solutions 
were calculated much more recently by Galfi and Racz . 
Their analytical solutions were obtained for a reaction of 
second order in an infinite domain with stepwise initial 
conditions, in the asymptotic limit of long times. Under 
these conditions, exact scaling laws were derived show- 
ing, for instance, that the position of the maximum of 
the production rate grows as t 1 / 2 while the width of the 
reaction production zone grows as i 1 / 6 . Furthermore, the 
approach to the asymptotic state was studied by Taitel- 
baum et al. and they observed a slow approach to the 
asymptotic state, which was modeled by using perturba- 
tion analysis and could not be described by the long time 
asymptotics. 

At the experimental level, one is often confronted with 
the additional presence of advection in the fluid. For 
this reason, reactions are often studied in small capillary 
tubes where no large scale flow is present. Such sys- 
tems have been used 0, Q to study a reaction process 
with initially separated reactants, validating the power 
law predictions of Galfi and Racz at asymptotically long 
times. Furthermore, recent experiments in porous me- 
dia have shown very good agreement between theory 
and experiment in the asymptotic regime, in the case of 
one diffusing reactant into a solid substrate. 

Here, we address the more general problem of two dif- 
fusing, initially separated, reactants in a regime that 
spans both the initial transients and the asymptotic 
state. In contrast with previously published work 0,13, 
we turn our attention to the reaction product and its con- 
centration, instead of the production rate. This analysis 
takes into account the diffusion of the product as well as 
the reagents. We compare numerical and experimental 
measures of the reaction region, beyond the scaling of the 
different quantities, to study the shape of the concentra- 
tion profiles. 

Microfluidic circuits are used to achieve low Reynolds 



numbers (Re ~ 1), thus removing any effects of advec- 
tion without resorting to gels or other media. We believe 
that the detailed comparison of the shape of the con- 
centration profiles between simulations and experiments 
allows one to assess to what extent, from the dynam- 
ical viewpoint, modeling a reaction by a kinetics of a 
prescribed order is accurate. This assessment is useful 
in view of modeling more complex behaviors in laminar, 
chaotic, or turbulent flows. In this regard, microsystems 
offer the opportunity to accurately compare theory and 
experiment since simple geometries can be achieved and, 
owing to the low Reynolds number at hand, the hydrody- 
namics is well controlled. Diffusion is the only transport 
mechanism leading the reagents to react together. 

The channel geometry we consider here precisely allows 
the study of a reaction-diffusion process with stepwise ini- 
tial conditions. The micro-reactor is T-shaped, similar to 
that used in Ref. |6(. In contrast with these authors, we 
use shallow channels which favor two-dimensionality as 
discussed below. We will further infer from the present 
work a practical method for measuring the kinetics of 
chemical reactions in the sub-millisecond range, a regime 
that is out of the reach of traditional stopped-flow tech- 
niques 0. 

Theoretical background: Let us consider a second or- 
der irreversible reaction of the type A + B — ► C, taking 
place in a two-dimensional space, within a flow driven at 
uniform velocity U, as shown in Fig. ^ 
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FIG. 1: Geometry of the theoretical problem. 



In the stationary state, after neglecting the derivative 
along X in the Laplacian operator, the reaction-diffusion 
equations read : 
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where A and B are the concentrations of the two reac- 
tants A and B, D a and Db are the corresponding diffusion 
coefficients, and k is the chemical rate constant. These 
equations are identical to the typical reaction-diffusion 
equations, if we replace variable X by Ut, where t is the 
time. This system must be supplemented by an equation 
governing the evolution of product C, given by : 
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where C is the concentration of product C, assumed in- 
ert, and D c is the corresponding diffusion coefficient. As 
initial conditions, we take the reagents to be separated 
with constant densities for both Y < (A = Ao and 
B = 0) and Y > (A = and B = B ). We thus have 
here a six dimensional parameter space which consists of 
the three diffusion coefficients of the two reagents and the 
product, the initial concentration levels of the reagents, 
and the reaction rate coefficient. 

The system of Eqs. ^ can be non-dimensionalized by 
introducing a characteristic length scale £, similar to the 
one introduced in Ref. Q 
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This length is used to non-dimensionalize X and Y as: 
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where D = ^D a Djj. Using these new dimensionless vari- 
ables and introducing the dimensionless concentrations 

a = — — and o = -— , 
A Bo 

the governing equations for the problem become: 
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in which we define 
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This system must also be supplemented by an equa- 
tion prescribing the evolution of the product concentra- 
tion. Although the problem can be solved in the general 
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FIG. 2: geometry of the microreactor. 



case, we restrict ourselves to the particular case where the 
third diffusion coefficient, D Cl is equal to D a (assumed to 
be the smaller than -Dfc). The equation for the evolution 
of the product concentration thus reads, in dimensionless 
form: 



dc 
dx 
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with c = C I \J AqBq- Note that the reaction rate has 
disappeared in Eqs. 21 an< i EH This means that chem- 
ical kinetics does not play any dynamical role: chang- 
ing k rescales the variables of the problem but does not 
change the structure of the solution. A consequence is 
that the number of parameters is reduced to two, mean- 
ing that the full set of solutions can be represented in a 
two-dimensional parameter space. A knowledge of both 
parameters \ an d P 1S necessary to compare simulations 
and experiments, if one is interested in both the transient 
and asymptotic regimes. 

Description of the experiment: In our experiments, we 
used a T-shaped microreactor, whose geometry is rep- 
resented in Fig [3 This geometry is similar to the one 
introduced in Ref. 0- 

The channels were chemically etched in a glass wafer, 
which was then anodically bonded to a silicon wafer 
to create a rectangular cross-section. Throughout the 
experiments discussed here, the channel depth was 10 
fim and the width was 200 /im. The lengths of the T 
arms were I cm while the main channel was 2 cm long. 
Two reactant species are introduced through the two top 
branches of the T and they react as they flow side by side 
in the main channel. We study the reaction of Calcium 
ions as reactant A and Calcium Green (CaGreen) from 
Molecular Probes @ as reactant B (as used, e.g. in 0). 
Ca 2+ ions bind with the tracer Calcium Green, signif- 
icantly increasing its fluorescence |lCj . Epifluorescence 
microscopy is used to measure the fluorescence emission, 
yielding spatial images such as Fig. [3J The flows are 
driven by syringe pumps which also provide measure- 
ments of the flow-rates in each separate branch; those 
ranged between 30 and 100 nl/mn. Owing to the aspect 
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FIG. 3: Fluorescent microscopy image of the binding of 
Ca 2+ ions with Calcium-green marker in a T-channel. The 
channel width is 200 /im, the depth is 10 /im. Here, the ve- 
locity in the main channel is U = 0.083 cm/s. 



ratio of the micro-channel, the velocity profile is parabolic 
across the channel depth and is essentially uniform along 
Y direction. Furthermore, the low Reynolds number en- 
sures that the fluid flow is invariant in the x-direction 
after some entrance length on the order of the channel 
thickness. 

In the conditions in which we operated (i.e. low flow- 
rates and thin channels) the tracer concentration can be 
considered as uniform across the channel depth. As a 
consequence, the concentration fields can be treated as 
two-dimensional and Eqs.0]and[S]can be used to describe 
the system. We tested the two-dimensionality property 
by performing a series of diffusion experiments (with no 
reaction) of fluorescein in water. A fluorescein solution 
and pure water are driven from the two top arms re- 
spectively. In the main channel, they form a diffusing 
interface which progressively broadens downstream. In 
thick channels, the broadening was found to display the 
three-dimensional effects discussed in Refs. and [TT| . 
where the diffusion region scales as a; 1 / 3 . In the 10 /zm 
channels, 2D diffusion theory to agrees well with experi- 
ment, with a scaling of the diffusion width as x 1 ^ 2 . 

The comparison between the 2D theoretical model and 
the experiment is further justified by the following rea- 
soning: one may consider that the reacting zone is homo- 
geneous across the channel depth in the major part of the 
region we analyze, for sufficiently thin channels. This ar- 
gument relies on estimations of the diffusive time across 
the channel thickness and on the measurements with flu- 
orescein discussed above. We are thus working mostly in 
the so-called Taylor- Aris regime where the parabolic ve- 
locity profile enhances the apparent diffusion coefficient 
in the flow direction One may show however that 

the streamwise diffusion term is still negligible compared 
to the transverse term, since the streamwise gradients 
are very weak. These two features justify the use of the 
model of Eqs.0]to compare with our experiment. 

The diffusion coefficient of Calcium ions is equal to 
D Ca2+ = 7.9 x 10~ 6 cm 2 /s 0. The diffusion coefficient 
of the CaGreen was unknown so we measured it using the 
method described above, replacing the fluorescein solu- 
tion with a Calcium Green solution: We checked that the 



width of the diffused zone increases as the downstream 
distance raised to the power 1/2, as expected from two- 
dimensional theory. By fitting the measurements of the 
width, we obtained an estimate for the Calcium Green 
diffusion coefficient equal to 3.2±0.3 x 1CP 6 cm 2 /s. 

Throughout the experiments, we kept the concentra- 
tion of Ca 2+ constant at [Ca 2+ ]=1 mM, while the con- 
centration of CaGreen was varied in the range 0.1 < 
[CaGreen] < 0.5 mM. Therefore, the value of (3 varied 
between 1.4 and 3.2, while x was fixed at \ = 1.57. 

Comparison between theory and experiment: A typical 
concentration field is shown in Fig EI The black region 
corresponds to the Calcium solution, the grey one to the 
CaGreen solution, and the brighter region in the center 
to the bound Ca-CaGreen complex. The background flu- 
orescence, visible in the lower part of the channel, is due 
to the fluorescence of unbound CaGreen. The asymmetry 
in the reaction zone is due to a different rate of diffusion 
of the reagents, because of their different initial concen- 
trations and different diffusion coefficients. This can be 
physically explained by noting that the species with the 
higher diffusion coefficient will diffuse further; further- 
more the species with the larger initial concentration will 
diffuse faster due to the large gradient it experiences. In 
other words, (3 and \ both determine the shape of the 
reaction zone. 

We were not able to obtain reliable information on the 
kinetic coefficient k of the Ca 2+ -Calcium Green reaction. 
We thus took k as a free fitting parameter, and adjusted 
it so as to optimize the agreement between theory and 
experiment. The best fits were obtained for k = 1.0 ± 
0.47 x 10 6 dm 3 /(mol s), or a length scale I = 0.84 /jm for 
f3 = \/2 and I = 1.0 ^m for (3 = 2. These values of k 
and £ are used throughout this letter in comparing the 
experiments and numerics. 

Figure 0] represents the concentration profiles of the 
Ca-CaGreen complex (dots) measured at various loca- 
tions downstream, for (3 = \J~2 and 2. The two locations 
we consider here are X = 278 and 1128 /im (for f3 = \/2 
these correspond to x = 32 and 216; for (3 — 2, x = 23 
and 153). The plots we show are representative of the 
evolution of the concentration profiles along the down- 
stream axis. The profiles are asymmetric with respect 
to Y = 0; in particular the maximum of the concentra- 
tion is not located at the channel center but is displaced 
towards the side of the CaGreen solution. For increas- 
ing X, the reaction zone broadens. The solid lines on 
Fig. 0] represent the value of c obtained by numerical so- 
lutions of Eqs.^]and|Sl by a standard one dimensional 
finite-difference scheme with time stepping on a personal 
computer, using the values of the parameters previously 
established. The offset on the right hand side of the ex- 
perimental curves is due to the background fluorescence 
of the unbound CaGreen (reactant a). We match this 
offset in the numerical curves by plotting ha + c rather 
than just c, where h a normalization constant that gives 
the correct offset. This does not affect the fitting param- 
eters but matches the experimentally measured quantity, 
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i.e. the total fluorescence. We observe good agreement 
between theory and experiment. 
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FIG. 4: Experimental measurements of the fluorescent in- 
tensity (dots) compared with numerical solutions of c (solid 
lines) at different downstream locations for f3 — v2 and 2. 
U — 0.334 cm/s for both cases. 

A quantitative measure of the agreement between the 
simulations and experiments is obtained by looking at 
some representative quantities for the profiles. In par- 
ticular, we compare the downstream evolution of three 
quantities: The position of the maximum, the width 
at half- height, and the value of the maximum (Fig. [SJ. 
These evolutions do not follow simple power laws in the 
full range of x we consider. However, for x > 30 dimen- 
sionless units, the power laws are better defined. The 
useful scaling range in the experiments is further limited 
(at X ~ 1200/im) by lateral wall effects, thus limiting the 
range over which a power law can be observed; although 
the scaling range is too small to obtain a clear exponent, 
our measurements are consistent with the long time limit 
analytical theory of Galfi and Racz pj. The maximum 
concentration shows a saturation effect at large distances, 
i.e above x ~ 300. These measurements can be compared 
with the numerics, plotted as full lines. Fig[5]shows there 
is good agreement between the theory and the experi- 
ment. A similar agreement is obtained different values of 
(3 and flowrates U we considered. 



In conclusion, we demonstrate that two-dimensional 
reaction-diffusion equations with stepwise initial condi- 
tions define not only a theoretical problem, but also a 
physical situation which can be achieved experimentally. 
Microfluidic reactors are particularly convenient in order 
to produce such physical situations. Conversely, 2D reac- 
tion diffusion equations can be used to determine, from 



experimental measurements, the relative concentration of 
the reactants (micro-titration), the diffusion coefficients 
or the chemical kinetics. In the present case, we demon- 
strate a determination of the kinetic coefficient k in order 
to match measured physical quantities (position, width 
and concentration of the reactant) with numerical predic- 
tions. The value of k thus obtained for our system cor- 
responds to a reaction half-time of 1 ms, which is faster 
than the time scales that can be resolved in a stopped- 
flow apparatus. This indicates that the present approach 
can be used for analyzing fast chemical kinetics. 

The authors acknowledge technical help from Bertrand 
Lambollez, Jean Rossier and Thibault Colin, and useful 
discussions with Martin Bazant. 
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FIG. 5: Evolution of the position of the maximum concentra- 
tion, width, and value of the maximum with X. The + signs 
are the experiments, while the solid lines are the numerical 
simulations. Here, f3 = V2, £ = 0.84 /im, U = 0.334 cm/s. 
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